Dynamic structure factor of a superfluid Fermi gas 
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We describe the excitation spectrum of a two-component neutral Fermi gas in the superfluid 
phase at finite temperature by deriving a suitable Random-Phase approximation with the technique 
of functional derivatives. The obtained spectrum for the homogeneous gas at small wavevectors 
contains the Bogoliubov- Anderson phonon and is essentially different from the spectrum predicted 
by the static Bogoliubov theory, which instead shows an unphysically large response. We adapt 
the results for the homogeneous system to obtain the dynamic structure factor of a harmonically 
confined superfluid and we identify in the spectrum a unique feature of the superfluid phase. 
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I. INTRODUCTION 



The techniques of atom trapping and cooling -which have led to the realization of Bose-Einstein condensation in alkali 
gases are currently being employed to cool also the fermionic isotopes ^"^K and ^Li Dilute gases of fermionic 
atoms -with attractive inter-particle interactions are predicted to undergo a superfluid transition at lo-w temperatures 
(Tsup ^ Tp, -where Tp is the Fermi temperature). The simplest mechanism envisaged is a s--wave pairing |^,^ -which 
can be obtained, compatibly -with the Pauli principle, bet-ween atoms belonging to f-wo different internal states. The 
realization of a t-wo-component Fermi gas in the superfluid state may provide a ne-w physical system to study. Its 
properties are expected to be different from those of superfluid '^He, -which has a p--wave pairing and is not in the 
dilute regime, and from conventional charged superconductors -which have an excitation spectrum dominated by the 
Coulomb interaction Q and only weakly modified by the superfluid transition. 

An important issue for future experiments is to identify a clear signature of the superfluid transition in an atomic 
Fermi gas. Contrary to the case of atomic Bose-Einstein condensates, for fermions the superfluid transition affects 
only slightly the density proflle and the internal energy of the gas . A first idea is to measure the pair distribution 
function of the atoms, e.g. by using a laser probe beam j7|. A second idea is to look at the dynamical properties, 
-which are expected to be dramatically modified by the transition. Several proposals have been put for-ward in this 
direction, such as Cooper-pair breaking via a Raman transition Q , measurement of the moment of inertia of the cloud 
1^ , and excitation of collective modes in a harmonic trap by modulation of the trap frequencies [p^ or rotation of the 
axis of the trap . 

In this -work we suggest to identify the superfluid phase through the measurement of the bulk excitations of the gas, 
i.e. excitations -with a -wavelength smaller than the spatial extension of the atomic cloud. This is complementary to the 
proposals in as it deals -with high energy excitations in a quasi-homogeneous system at arbitrary temperature, 

and is not irrealistic from the experimental point of vie-w, since efficient Bragg scattering techniques have already 
been successfully used to measure the excitation spectra of Bose condensates lia] . 

We obtain the excitation spectrum of the fluid in the dilute regime by employing the Random-Phase (or Time- 
Dependent Hartree-Fock-Gorkov) Approximation (RPA), which we derive explicitly for the two-component system 
of present interest by the technique of functional derivatives ^3 14 . We do not use here the usual static Bogolubov 
approximation [ p^ for two main reasons: (i) physically its excitation spectrum has a gap and therefore ignores the 
branch of phonon-like excitations (Bogolubov- Anderson phonon, [[l6| - |l8[ ) expected on very general grounds to show 
up in homogeneous neutral superconductors with short range interactions and (ii) the density-density response 
function obtained within this approximation shows an unphysically large response at small wavevectors and fails to 
satisfy the /-sum rule, which is a requirement deriving from the local particle conservation law. 

We find that the RPA spectrum of a homogeneous two-component Fermi gas in the superfluid phase possesses the 
continuum of particle-hole excitations and a peak corresponding to the Bogoliubov- Anderson phonon; it satisfies the 
/-sum rule and includes naturally the Landau damping of the phonon due to the interplay with thermal excitations. 
We adapt the results obtained for the homogeneous system to describe harmonically confined gases by means of a 
local-density approximation, which predicts a broadening of the spectrum by taking into account the inhomogeneity 
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of the density profile. We predict that the Bogoliubov- Anderson phonon, which is the main feature of the spectrum 
in the superfLuid phase, would remain visible even in the trapped cloud. 

The structure of this paper is as follows. In Sec. II we derive the density response function of the fluid in the RPA. 
In Sec. Ill we obtain the spectrum of density fluctuations first in the homogeneous case and then for an harmonically 
trapped gas. Finally, Sec. IV gives a summary of our results and offers some concluding remarks. 

II. RANDOM-PHASE APPROXIMATION 

We describe a two-component atomic Fermi gas by the following Hamiltonian: 

+ ^ E /'^'^ / dV'v31(r)V'^(r')z.„Mr-r')^;3(r')V'Jr) . (1) 

{a,/3} = T,i 

The fermionic field operators V'a(r) satisfy the usual anticommutation relations; the interactions are considered only 
in s-wave between fermions in different internal states and arc modeled by the inter-particle potential Va_fi(y\ — r2) — 

,^ ^ X ^Til*"! r2) Y -vvith W|j^(ri — v-2) — v^^&ivi — r2)5ri2(''i2')- This model interaction potential, known 

as the Fermi pseudopotential, leads to a divergence-free BCS theory For the sake of generality we have included 
the presence of the external confinement via the trapping potentials K!xt,a(r) and the possibility of having different 
numbers of atoms in the two components; however in the following we shall restrict to the derivation of the equations 
in the homogeneous system and in the symmetric case = = N/2, which is the most favourable for the formation 
of Cooper pairs Q . The effect of trapping potential present in a realistic experiment will be included later on with a 
local density approximation. The external perturbing field which is necessary to generate the total density response 
has been introduced as U (r, t). Physically U represents the action of the probe applied in a real experiment; it may be 
a time-dependent perturbation applied to the magnetic trap |Q , a probe laser beam or a test particle |^ . Here 
we have assumed that the same potential U acts in the same way on both components, and we shall determine the 
perturbation Sp{r,t) on the total density induced by the probe potential U{r,t), assuming that the gas is initially at 
thermal equilibrium with a temperature T. We restrict to the linear response regime, where the density perturbation 
is a linear functional of the probe potential U expressed through the density-density response function x(ri, ii, r2j ^2)1 
a function of two position vectors ri.2 and of two time variables ti_2- 

6p{YiM)^ J d^r2dt2 xiri,ti,r2,t2)U{r2,t2). (2) 

In this section we explain how to calculate this response function in the Hartree-Fock-Gorkov approximation. We 
obtain general equations valid for an arbitrary trapping potential 14xt(r)- We then solve these equations explicitly 
for a spatially homogeneous gas at thermal equilibrium, where x is a function of ri — r2 and ti — t2 only. 

To proceed with the derivation of the response function x, we follow the imaginary time Green's function technique 
of the book of Kadanoff and Baym jl^. One first defines the two by two matrix of (normal and anomalous) Green's 
functions in imaginary times 

^(^'2) =(§;[;; 2) g;;(;;2))--™i)*^(2)) (3) 

where T is the time-ordering operator, "^(l) — i ^t^i! )' 'I'U^) = ("01(2), ^1(2)1, (...) indicates the average over 



the state of the system in the presence of the perturbing field U and (1, 2) stands for (ri, iri , r2 , iT2), where ti and T2 



are real quantities. For more details on the imaginary time technique, we refer to reference 13 1. We simply note that 
the various functions considered here can be obtained for real times by analytic continuation of their imaginary time 
values. From the equation of motion for the field operator in imaginary times one derives the generalized Dyson 
equation for G: 

G(l, 2) = Go(l, 2)+ [ d3 [ M Go(l, 3)S(3, 4)G(4, 2) + [ d3 Go(l, 3)M^(3)G(3, 2) , (4) 
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where the 2x2 matrix Go(l,2) is the solution of the equations of motion in absence of the interactions, W{1) 
fU{l) 



-{/(I) 



is the 2x2 matrix of external field and 2) is the 2x2 matrix of self-energies. Since we want to 



describe a dilute system, we work in the mean-field Hartree-Fock-Gorkov symmetry breaking approximation, where 
the self-energy reads 

V n 9^ An 9^ / (^l(l)V'i(l)) \^ no^f « GTi(l,2)\ 

i]Hc(i,2) = 5(i,2).n(^ -(v3t(i)^,(i)> ' Hg,Ki,2) )■ 

Here we have used the notations (5(1,2) = 5{ri — f2)S{Ti — T2)/i and Ufj(l,2) = U|j^(ri, r2)(5(Ti — T2)/i. We remark 
that the Fock contribution is zero since the interaction takes place only between particles with opposite spins and 
= in the considered state of the system, hence the vanishing diagonal in the last term of Eq. 
Following a standard approach [ ^"slJlill , we then obtain the density response matrix in RPA by taking the functional 
derivative of the Green's function with respect to the external field U: we define the generalized response matrix 

as the two by two matrix i(l,2,3) — CT3(5G(1, 2)/(5J7(3), where cts is the third Pauli matrix ^ q ^1 ) "^^"^ matrix 

giving the physical response is obtained from the limit L{1, 2) = £(1, 1^, 2). The density-density response function x 
is simply the trace over the two spin components of the response matrix: 

X(l, 2) - Tr 2) - L^^{1, 2) + Lu(l, 2). (6) 

The equation for the density response in the Random-Phase Approximation is obtained by the functional derivative 
of the Dyson equation, Eq. (P, where the approximation for the self energy has been employed. This yields: 

L(1,2) = L"(1,2)-|-^ f d3L°(l,3)x(3,2)- f d3 f d4 G(l, 3)M(3, 4, 2)G(4, l)«u(3, 4) . (7) 



Physically i*'(l,2) is the response matrix of the gas in the static Bogoliubov approxi- 



Here we have introduced the 2 x 2 matrices G(l,2) = cr3G(l,2), L°(l,2) = G(l, 2)G(2, 1), and M(l,2,3) 
Ln{l,2,3) 
iiT(l,2,3) 

mation, so that the reference system of the RPA is not the ideal gas but the Bogoliubov gas of quasiparticles. 

It is possible to display the diagrammatic structure of Eq. by separating out the "proper" part L{1,2) of the 
density response. We have therefore an equation which sums the bubble diagrams, 

L{1, 2) = Z(l, 2) + ^v^^ f d3 L{1, 3)x(3, 2) , (8) 
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and an equation which defines the bubble as a sum of all the ladder diagrams 



Z(1,2) = L"(1,2)- / d3 / d4G(l,3)«u(3,4)M(3,4,2)G(4,l) . (9) 



Here M = ( 

We now specialize the previous equations to the case of a spatially homogeneous gas and to the dilute limit A <C eF, 
where A is the gap and is the Fermi energy. All the response matrices depend then only on the relative spatial 
coordinates r = ri — r2 and the relative time coordinate t — ti — t2 and we introduce their double Fourier transforms 
with respect to r and t, e.g. 

L{q,iu)= y"dV(itL(ri=r,ii=t,r2 = 0,i2 = 0)e^(i'-"*). (10) 

The following equation is obtained for L(q, lj) from the solution of Eq. (||) together with the regularization of the 
contact potential: 

L(q,c.) = Aici,u)+vn-^^^^^^, , (11) 

where j4(q, w), i?rcg(q, i-^) and G(q, to) are complex functions of the frequency to be evaluated numerically. We remark 
that the assumption A <C ef has considerably simplified the treatment by allowing to introduce only three basic 
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functions {A, B and C) in place of six required by the exact treatment 
R stands for A, B 01 C and r stands for a, b, c, is given by 



The general expression i?(q, w), where 



(fik 



rL(k,q)(/(i?4 



(27r)3 

+rB(k,q)(l-/(£;^ 



f{E-)) 
fiE-)) 



1 



huj + E+ - + ir] 
1 



hio- (E+ -E^] 



IT] 
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fiLU - {E+ + E-) + h] 

where the "Landau" contributions ri(k, q) are given by 

ai(k,q) = (1 + (e+e- - ^^)/E+E^)lA 
6L(k, q) = (1 - (e+C- + ^^)/E+E^)lA 
CL(k,q) = A(l/i;+-l/S_)/8 

and where the "Beliaev" contributions rB(k, q) are given by 

as(k,q) = (1 - (e+e- - A2)/i?+i?_)/4 
6B(k,q) = (1 + (e+e- + A2)/i?+i?_)/4 
cs(k,q) =-A(l/i?+ + !/£;_ )/8 



a 



hw + E+ + E- + iri 



(12) 



(13) 
(14) 
(15) 



(16) 
(17) 
(18) 



for A(q, w), B{q,uj) and C(q, w) respectively. The function i?(q, w) as defined in Eq. (|T^) presents an ultraviolet 
divergence originating from the choice of a contact interaction potential. This divergence is removed in a systematic 
way by the use of the pseudopotential, which amounts here simply to subtracting the most diverging contribution in 
the form l/(2^k): 



/d^k /I 



(19) 



where V represents the principal value |2|]. We have used the notations ^± — fi,^(k± q/2)^/2m + U||n/2 — /i where n 
is the total equilibrium density of particles, E± = (^^ + A^)^/^ and f{E) — l/(exp(/3i?) + 1) is the Fermi distribution 
function at temperature T with f3 = l/{kBT). The parameter a equals 1 for A(q, cj) and -B(q, lj), while a equals — 1 
for C{(\,uj), and 77 is a positive infinitesimal. 

The final expression for the double Fourier transform of the density-density response function in RPA is obtained 
from Eq. (§) as: 



2£(q,^) 
1 - v^lLici.uj) 



(20) 



This can be contrasted with the density-density response function in the static Bogoliubov approximation, which 
leads to x(q,w) = 2A(q,w) 



III. THE SPECTRUM OF DENSITY FLUCTUATIONS 



A. Homogeneous system 

Before displaying the fully numerical solution of Eqs. ( pi] ) and ([20|), we analyze some limiting cases. At temperatures 
higher than the BCS transition temperature we have i?(q, lo) — 0, C(q, lo) — Q and A(q, lo) = Xo(qj where Xo(qi ^) 
is the well-known Lindhard function for the response of a ideal Fermi gas (see for example [^). The RPA equation ( pl| ) 
reduces to the usual expression x = 2xo/(l — ''^TiXo), the factor two being due to the two spin components of the gas. 
In the case of repulsive interactions {i.e. v^i > 0) the equation shows a pole corresponding to the zero sound, while 
no well-defined collective excitation is stable in the case of attractive interactions. 

At zero temperature, in the limit of small q and m it is possible to estimate analytically the expression for the 
density response function; to lowest order in q wc obtain 

[uj + iviy - c%q\l + 2kFa^i/TT) 
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where cb = vp/VS is the sound velocity predicted by BogoUubov, Af{eF) = rnkp/n'^h'^ is the density of states at the 
Fermi level and an is the scattering length such that = Airh'^an/m. The pole yields a phonon-like excitation, 
corresponding to the Bogoliubov- Anderson sound for this system. Eq. ( pl| ) holds approximately until the phonon is 
stable, that is before it meets the continuum of quasiparticle-quasihole excitations, which has a threshold energy of 
2A. Evidently the RPA expression is valid in the dilute limit fci?|a||| ^ 1. 

It is easily checked that in the long-wavelenght limit the Bogoliubov-Anderson sound exhausts the /-sum rule 
— J CL;Imx(q, cj) = nnq'^/m, where n is the total density of the gas. A more general proof can be obtained 
by noticing that the RPA, being equivalent to the time-dependent Hartree-Fock-Bogoliubov theory, automatically 
satisfies the continuity equation and hence the /-sum rule. On the contrary, the static Bogoliubov approximation 
results to be bad in the limit q ^ 0: from Eq. ( [12| ) we estimate that A(q ^ 0, w) oc oj"'^/^ in the high frequency limit, 
yielding an infinite contribution to the first moment integral. 

We now turn to the presentation of numerical results. Rather than plotting the complex quantity t^) we have 
chosen to represent the spectrum of total density fluctuations of wavevector q, given by the dynamic structure factor 
S'(q, uj). On an experimental point of view the dynamic structure factor can be accessed via the rate of the scattering 
events of a probe particle by the gas leading to a momentum exchange hq and to an energy exchange hco between 
the probe particle and the gas; on a theoretical point of view the dynamic structure factor of the gas is related to 
Imx(q, ti^) by the fluctuation-dissipation theorem |^ : 

S{q, Lo) = -(2ft/n)(l - exp(-/3fii^))-ilmx(q, to). (22) 

Figure ^ shows the spectrum of a homogeneous superfluid at zero temperature as resulting from the full RPA calcula- 
tion, compared to the predictions of the static Bogoliubov approximation: it is evident that for small q {hq < 2A/cb) 
the Bogoliubov approximation yields an unphysically large response (Fig. 1 (a)). For larger q the Bogoliubov-Anderson 
phonon falls in the continuum of quasi-particle quasi-hole excitations, and the two approximations yield almost the 
same result (Fig. 1 (b)), which is also close to the ideal-gas solution. 



B. Harmonically trapped system 

We turn now to the situation where the particles are subject to an external harmonic conflncmcnt. We assume that 
the conflning potential has the same action on both spin components, this is indeed the case in a laser induced trap. 
For simplicity we further assume that the resulting trapping potential is isotropic so that 

Vc^t,^{r) = 14xta(^) = ^rnflW (23) 

We consider bulk excitations of the harmonically confined cloud induced by a probe potential U{r, t) of wavevector q 
and frequency uj. We characterize the density response of the gas by the dynamic structure factor calculated in the 
local-density approximation: 

SLDAin, u;) = ^ _ J d\ Imx(q, ^x{r), A(r)) , (24) 

where x is the density response function derived in the previous section for the homogeneous system. This local- 
density approximation is valid for q > l/R, where R is the radius of the cloud, since it does not take into account 



surface modes |10 and in-gap single-particle excitations pq |. The same approach has already described successfully 
an experiment on trapped Bose-Einstein condensates, where the experiment has directly measured the function 
Imx(q,c^) (X [S{ci,Lj) ~ 5(-q, -co)] 0]. 



The position-dependent chemical potential and gap in Eq. (24) are defined as fj,(r) = ^ — T4xt('') and A(r) = A[n(r)]. 
The chemical potential /i of each spin component is determined by the normalization condition N = J d^rn{r) where 
N is the total number of particles in the gas. The equilibrium total density profile n{r) and the gap A(r) are obtained 
first by numerical solution of the BCS equations in the homogeneous system 

^{Kiv(£^k) + kkni-/(i?k)]} , (25) 

(2^)3 \ 2E^ UeJ/~ «U ' ^ ' 
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and then by employing the Thomas- Fermi approximation (TFA) to take into account the inhomogeneity due to 
the external confinement. 

The main effect of the external confinement is a broadening of the spectrum, which is due to the inhomogeneous 
distribution of the density in the trap. This is illustrated already by a simple analytic expression for the dynamic 
structure factor at zero temperature in the low q limit: integrating the imaginary part of Eq.(pl|) over space yields to 
lowest order in /ci?|a|J 



Here, v4g[7i] = rnkpCBq/'^'^fi, = csg and we have adopted the rescaled units uj = oj/ft, q — hq/\/mhn and 
p, — fi/hQ. For the density profile we have taken the TFA expression n{r) — (2m{ii — 14xt('^))/ft^'^^^/37r^, which 
neglects the Hartree mean field effect but turns out to be a good approximation in the dilute limit B. 

The numerical results at finite temperature and wavevector, presented in Fig. ^ show in the low temperature 
spectrum a peak corresponding to the Bogoliubov-Anderson phonon and include a high-frequency tail due to the 
contribution of multi-particle excitations. With increasing temperature, the asymmetric feature due to the Bogoliubov- 
Anderson phonon becomes less marked and disappears at T ~ A/fc^, when quasi-particle quasi- hole pairs are easily 
excited by thermal fluctuations. 



IV. CONCLUSION 



In this paper, we have derived for a two-component spin-polarized Fermi gas a generalized Random-Phase Approx- 
imation to describe the excitation spectrum in the superfluid phase. We have shown that, contrary to the case of 
bosonic systems, the predictions of this theory - valid in the limit A;_F|af J ^ 1 ~ are essentially different from those 
of the "static" Bogoliubov theory, which instead yields an unphysically large signal at small wavevectors, due to the 
lack of local particle conservation. 

The possible experiments that have motivated this theoretical work are light scattering or scattering of test 
particles ]28| , pT[| by a two-component Fermi gas stored in a dipole trap. The outcome of this type of experiments is 
described by the dynamic structure factor. We have therefore employed the results of the homogeneous system to 
predict in a local density approximation the dynamic structure factor of a harmonically trapped superfluid at finite 
temperature, and we have shown that the Bogoliubov-Anderson phonon - main feature of the superfluid phase in 
the spectrum of the homogeneous gas - would appear also in the response of the trapped system, as an asymmetric 
peak. We have investigated the effects of the temperature on the shape of the response, showing that the Bogoliubov- 
Anderson phonon should remain visible up to a temperature T ~ A/ks- The observation of the Bogoliubov-Anderson 
phonon in the response of the gas to a probe beam may therefore provide a way to detect the presence of the superfluid 
phase in the experiments on alkali Fermi gases. 

Our general RPA equations are also suitable for a full description of the inhomogeneous system without local density 
approximation, thus allowing in principle to take into account the discrete nature of the eigenmodes of the trapped 
gas. This would complete the static Bogoliubov treatment already performed in |pj|. 
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FIG. 1. Dynamic structure factor of a homogeneous superfluid Fermi gas at T = for q = 0.054 kp (a), and for q — 0.54 
kp (b), as predicted by the RPA (sohd line) and by the static Bogoliubov approximation (dashed hue). The results for a 
non-interacting Fermi gas (long-dashed line) are also shown, kp is the Fermi wavevector and ep = h^kp/2m is the Fermi 
energy of the non-interacting Fermi gas at the same density. The parameters chosen are A = 0.065eF and kpa-^i = —0.04. 
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FIG. 2. Dynamic structure factor of a confined superfluid Fermi gas for hq — 4.15 Pho at different temperatures; T — 
(solid line), kgT = 0.014 Ep (short-dashed line), hgT = 0.028 Ep (long-dashed line) and kgT = 0.042 Ep (dot-dashed 
line). The critical temperature for the BCS transition is fc^T = 0.047 Ep. The analytic expression for the broadening of the 
Bogoliubov-Anderson phonon in the low q and zero temperature limit, Eq. (^^, is given in dotted line. The inset shows the 
gap function in units of Ep as a function of the radial coordinate r/uho, as obtained from the local-density approximation. The 
different line-styles correspond to the different temperatures of the main figure, pho = VTimQ,, aho = \/h/mQ, and Ep is the 
Fermi energy of a harmonically trapped non-interacting Fermi gas with the same number of particles. The parameters used 
are N — 8 x 10^, f2 = 27r x 27.2 s~^ and = —2160 ao where ao is the Bohr radius. 
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